1 What is this part of the hands-on about?

Questions:

Data analysis in R using the tidyverse meta-package

Objectives:

Describe the purpose of the dplyr and tidyr packages.
Describe several of their functions that are extremely useful to manipulate data.
Describe the concept of a wide and a long table format, and see how to reshape a data frame from one format to the other one.

Keypoints:

Tabular data in R using the tidyverse meta-package

2 Data manipulation using `dplyr` and `tidyr`

Bracket subsetting is handy, but it can be cumbersome and difficult to read, especially for complicated operations.

Some packages can greatly facilitate our task when we manipulate data. Packages in R are basically sets of additional functions that let you do more stuff. The functions we’ve been using so far, like str() or data.frame(), come built into R; Loading packages can give you access to other specific functions. Before you use a package for the first time you need to install it on your machine, and then you should import it in every subsequent R session when you need it.

The package dplyr provides powerful tools for data manipulation tasks. It is built to work directly with data frames, with many manipulation tasks optimised.
As we will see latter on, sometimes we want a data frame to be reshaped to be able to do some specific analyses or for visualisation. The package tidyr addresses this common problem of reshaping data and provides tools for manipulating data in a tidy way.

To learn more about dplyr and tidyr after the workshop, you may want to check out this handy data transformation with dplyr cheatsheet and this one about tidyr.

The tidyverse package is an “umbrella-package” that installs several useful packages for data analysis which work well together, such as tidyr, dplyr, ggplot2, tibble, etc. These packages help us to work and interact with the data. They allow us to do many things with your data, such as subsetting, transforming, visualising, etc.

To install and load the tidyverse package type:

BiocManager::install("tidyverse")

# load the tidyverse packages, incl. dplyr
library("tidyverse")

3 Loading data with tidyverse

rna <- read_csv("../data/rnaseq.csv")

# view the data
rna

Notice that the class of the data is now referred to as a “tibble”.

Tibbles tweak some of the behaviors of the data frame objects we introduced in the previously. The data structure is very similar to a data frame. For our purposes the only differences are that:

It displays the data type of each column under its name. Note that <dbl> is a data type defined to hold numeric values with decimal points.
It only prints the first few rows of data and only as many columns as fit on one screen.

We are now going to learn some of the most common dplyr functions:

select(): subset columns
filter(): subset rows on conditions
mutate(): create new columns by using information from other columns
group_by() and summarise(): create summary statistics on grouped data
arrange(): sort results
count(): count discrete values

4 Selecting columns and filtering rows

To select columns of a data frame, use select(). The first argument to this function is the data frame (rna), and the subsequent arguments are the columns to keep.

select(rna, c(gene, sample, tissue, expression))

# select(rna, c("gene", "sample", "tissue", "expression"))

To select all columns except certain ones, put a “-” in front of the variable to exclude it.

select(rna, -tissue, -organism)

This will select all the variables in rna except tissue and organism.

To choose rows based on a specific criteria, use filter():

filter(rna, sex == "Male")

filter(rna, sex == "Male" & infection == "NonInfected")

Now let’s imagine we are interested in the human homologs of the mouse genes analysed in this dataset. This information can be found in the last column of the rna tibble, named hsapiens_homolog_associated_gene_name. To visualise it easily, we will create a new table containing just the 2 columns gene and hsapiens_homolog_associated_gene_name.

genes <- select(rna, gene, hsapiens_homolog_associated_gene_name)
genes

Some mouse genes have no human homologs. These can be retrieved using filter() and the is.na() function, that determines whether something is an NA.

filter(genes, is.na(hsapiens_homolog_associated_gene_name))

If we want to keep only mouse genes that have a human homolog, we can insert a “!” symbol that negates the result, so we’re asking for every row where hsapiens_homolog_associated_gene_name is not an NA.

filter(genes, !is.na(hsapiens_homolog_associated_gene_name))

5 Pipes

What if you want to select and filter at the same time? One way to accomplish this is using pipes, which are a recent addition to R. Pipes let you take the output of one function and send it directly to the next, which is useful when you need to do many things to the same dataset.

Pipes in R look like %>%.

In the above code, we use the pipe to send the rna dataset first through filter() to keep rows where sex is Male, then through select() to keep only the gene, sample, tissue, and expressioncolumns.

The pipe %>% takes the object on its left and passes it directly as the first argument to the function on its right, we don’t need to explicitly include the data frame as an argument to the filter() and select() functions any more.

rna %>%
  filter(sex == "Male") %>%
  select(gene, sample, tissue, expression)

Some may find it helpful to read the pipe like the word “then”. For instance, in the above example, we took the data frame rna, then we filtered for rows with sex == "Male", then we selected columns gene, sample, tissue, and expression.

The dplyr functions by themselves are somewhat simple, but by combining them into linear workflows with the pipe, we can accomplish more complex manipulations of data frames.

If we want to create a new object with this smaller version of the data, we can assign it a new name:

rna3 <- rna %>%
  filter(sex == "Male") %>%
  select(gene, sample, tissue, expression)

rna3

Challenge

Using pipes, subset the rna data to keep observations in female mice at time 0, where the gene has an expression higher than 50000, and retain only the columns gene, sample, time, expression and age.

6 Mutate

Frequently you’ll want to create new columns based on the values of existing columns, for example to do unit conversions, or to find the ratio of values in two columns. For this we’ll use mutate().

To create a new column of time in hours:

rna %>%
  mutate(time_hours = time * 24) %>%
  select(time, time_hours)

You can also create a second new column based on the first new column within the same call of mutate():

rna %>%
  mutate(time_hours = time * 24,
         time_mn = time_hours * 60) %>%
  select(time, time_hours, time_mn)

7 Split-apply-combine data analysis

Many data analysis tasks can be approached using the split-apply-combine paradigm: split the data into groups, apply some analysis to each group, and then combine the results. dplyr makes this very easy through the use of the group_by() function.

rna %>%
  group_by(gene)

The group_by() function doesn’t perform any data processing, it groups the data into subsets: in the example above, our initial tibble of 32428 observations is split into 1474 groups based on the gene variable.

We could similarly decide to group the tibble by the samples:

rna %>%
  group_by(sample)

Here our initial tibble of 32428 observations is split into 22 groups based on the sample variable.

Once the data has been grouped, subsequent operations will be applied on each group independently.

7.1 The `summarise()` function

group_by() is often used together with summarise(), which collapses each group into a single-row summary of that group.

group_by() takes as arguments the column names that contain the categorical variables for which you want to calculate the summary statistics. So to compute the mean expression by gene:

rna %>%
  group_by(gene) %>%
  summarise(mean_expression = mean(expression))

We could also want to calculate the mean expression levels of all genes in each sample:

rna %>%
  group_by(sample) %>%
  summarise(mean_expression = mean(expression))

But we can can also group by multiple columns:

rna %>%
  group_by(gene, infection, time) %>%
  summarise(mean_expression = mean(expression))

## `summarise()` has grouped output by 'gene', 'infection'. You can override using
## the `.groups` argument.

Once the data is grouped, you can also summarise multiple variables at the same time (and not necessarily on the same variable). For instance, we could add a column indicating the median expression by gene and by condition:

rna %>%
  group_by(gene, infection, time) %>%
  summarise(mean_expression = mean(expression),
            median_expression = median(expression))

## `summarise()` has grouped output by 'gene', 'infection'. You can override using
## the `.groups` argument.

Challenge

Calculate the mean expression level of gene “Dok3” by timepoints.

8 Reshaping data

In the rna tibble, the rows contain expression values (the unit) that are associated with a combination of 2 other variables: gene and sample.

All the other columns correspond to variables describing either the sample (organism, age, sex, …) or the gene (gene_biotype, ENTREZ_ID, product, …). The variables that don’t change with genes or with samples will have the same value in all the rows.

rna %>%
  arrange(gene) #sorts the dataframe

This structure is called a long-format, as one column contains all the values, and other column(s) list(s) the context of the value.

In certain cases, the long-format is not really “human-readable”, and another format, a wide-format is preferred, as a more compact way of representing the data. This is typically the case with gene expression values that scientists are used to look as matrices, were rows represent genes and columns represent samples.

In this format, it would therefore become straightforward to explore the relationship between the gene expression levels within, and between, the samples.

To convert the gene expression values from rna into a wide-format, we need to create a new table where the values of the sample column would become the names of column variables.

The key point here is that we are still following a tidy data structure, but we have reshaped the data according to the observations of interest: expression levels per gene instead of recording them per gene and per sample.

The opposite transformation would be to transform column names into values of a new variable.

We can do both these of transformations with two tidyr functions, pivot_longer() and pivot_wider() (see here for details).

8.1 Pivoting the data into a wider format

Let’s select the first 3 columns of rna and use pivot_wider() to transform the data into a wide-format.

rna_exp <- rna %>%
  select(gene, sample, expression)
rna_exp

pivot_wider takes three main arguments:

the data to be transformed;
the names_from : the column whose values will become new column names;
the values_from: the column whose values will fill the new columns.

Figure 1: Wide pivot of the rna data

rna_wide <- rna_exp %>%
  pivot_wider(names_from = sample,
              values_from = expression)
rna_wide

Note that by default, the pivot_wider() function will add NA for missing values.

Let’s imagine that for some reason, we had some missing expression values for some genes in certain samples. In the following fictive example, the gene Cyp2d22 has only one expression value, in GSM2545338 sample.

rna_with_missing_values

By default, the pivot_wider() function will add NA for missing values. This can be parameterised with the values_fill argument of the pivot_wider() function.

rna_with_missing_values %>%
  pivot_wider(names_from = sample,
              values_from = expression)

rna_with_missing_values %>%
  pivot_wider(names_from = sample,
              values_from = expression,
              values_fill = 0)

8.2 Pivoting data into a longer format

In the opposite situation we are using the column names and turning them into a pair of new variables. One variable represents the column names as values, and the other variable contains the values previously associated with the column names.

pivot_longer() takes four main arguments:

the data to be transformed;
the names_to: the new column name we wish to create and populate with the current column names;
the values_to: the new column name we wish to create and populate with current values;
the names of the columns to be used to populate the names_to and values_to variables (or to drop).

Figure 2: Long pivot of the rna data

To recreate rna_long from rna_wide we would create a key called sample and value called expression and use all columns except gene for the key variable. Here we drop gene column with a minus sign.

Notice how the new variable names are to be quoted here.

rna_long <- rna_wide %>%
    pivot_longer(names_to = "sample",
                 values_to = "expression",
                 -gene)
rna_long

Note that if we had missing values in the wide-format, the NA would be included in the new long format.

Remember our previous fictive tibble containing missing values:

rna_with_missing_values

wide_with_NA <- rna_with_missing_values %>%
  pivot_wider(names_from = sample,
              values_from = expression)
wide_with_NA

wide_with_NA %>%
    pivot_longer(names_to = "sample",
                 values_to = "expression",
                 -gene)

9 Exporting data

Now that you have learned how to use dplyr to extract information from or summarise your raw data, you may want to export these new data sets to share them with your collaborators or for archival.

Similar to the read_csv() function used for reading CSV files into R, there is a write_csv() function that generates CSV files from data frames.

Before using write_csv(), we are going to create a new folder, data_output, in our working directory that will store this generated dataset. We don’t want to write generated datasets in the same directory as our raw data. It’s good practice to keep them separate. The data folder should only contain the raw, unaltered data, and should be left alone to make sure we don’t delete or modify it. In contrast, our script will generate the contents of the data_output directory, so even if the files it contains are deleted, we can always re-generate them.

Let’s use write_csv() to save the rna_wide table that we have created previously.

write_csv(rna_wide, file = "data_output/rna_wide.csv")

Manipulating and analysing data with dplyr

Hands-on adapted from Introduction to R and Bioconductor

10 November 2022

1 What is this part of the hands-on about?

2 Data manipulation using `dplyr` and `tidyr`

3 Loading data with tidyverse

4 Selecting columns and filtering rows

5 Pipes

6 Mutate

7 Split-apply-combine data analysis

7.1 The `summarise()` function

8 Reshaping data

8.1 Pivoting the data into a wider format

8.2 Pivoting data into a longer format

9 Exporting data

Manipulating and analysing data with dplyr

Hands-on adapted from Introduction to R and Bioconductor

10 November 2022

1 What is this part of the hands-on about?

2 Data manipulation using dplyr and tidyr

3 Loading data with tidyverse

4 Selecting columns and filtering rows

5 Pipes

6 Mutate

7 Split-apply-combine data analysis

7.1 The summarise() function

8 Reshaping data

8.1 Pivoting the data into a wider format

8.2 Pivoting data into a longer format

9 Exporting data

2 Data manipulation using `dplyr` and `tidyr`

7.1 The `summarise()` function